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Abstract 

A first-principles-based technique is developed to investigate properties of Ba(Zr,Ti)03 relaxor 
ferroelectrics as a function of temperature. The use of this scheme provides answers to important, 
unresolved and/or controversial questions, such as: what do the different critical temperatures 
usually found in relaxors correspond to? Do polar nanoregions really exist in relaxors? If yes, do 
they only form inside chemically-ordered regions? Is it necessary that antiferroelectricity develops 
in order for the relaxor behavior to occur? Are random fields and random strains really the 
mechanisms responsible for relaxor behavior? If not, what are these mechanisms? These ab-initio- 
based calculations also leads to a deep microscopic insight into relaxors. 
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Relaxor ferroelectrics are characterized by some striking anomalous properties (see, e.g., 
Refs. [l-21] and references therein). For instance, they adopt a peak in their ac dielectric 



response-versus-temperature function whi 
cubic down to the lowest temperatures ll 



e they remain macroscopically paraelectric and 
. Furthermore, this dielectric response deviates 
from the "traditional" Curie- Weiss law [22| for temperatures lower than the so-called Burns 
temperature [2J. Other examples of anomalous properties include the plateau observed in 



their static, dc dielectric response at low temperature |23 
behavior 



24j, and the unusual temperature 



16] of the Edwards- Anderson parameter [25]. Determining the origin of these 



intriguing effects has been a challenge to scientists since the discovery of ferroelectric relaxors. 

The goal of this Letter is to report ab-initio-based calculations that not only reproduce all 
the aforementioned intriguing features but also offer a deep microscopic insight into relaxors. 

Practically, we decided to focus on a specific relaxor, namely disordered Ba(Zr .5Ti .5)O3 
(BZT) solid solutions (BZT is also fascinating because its parent compounds are rather dif- 
ferent: BaZrOs is paraelectric while BaTi03 is a typical ferroelectric). Here, we develop and 
use a first-principles-based effective Hamiltonian approach for which a detailed description 
is given in the Supplemental Material. The total energy of this effective Hamiltonian is used 
in Monte-Carlo (MC) simulations to compute finite-temperature properties of BZT alloys. 
We use 12 x 12 x 12 (8640 atoms) or 16 x 16 x 16 (20480 atoms) supercells in which the <Tj 
variables are randomly placed and kept fixed during the MC simulations, in order to mimic 
disordered BZT solid solutions. These two supercells provide similar results, which attest 
the convergency of the simulations. The temperature T is decreased in small steps from 
high temperature, and up to 10 6 MC sweeps are used to get converged statistical properties. 

Here, the Uj local soft-mode vectors in each 5-atom cell % (uj is directly proportional to 
the local electric dipole moment centered in cell i) and the homogeneous strain tensor tjh 
arising from the MC simulations indicate that Ba(Zr .5Ti .5)O3 bulk remains macroscop- 
ically cubic and non-polar for any temperature down to the lowest one investigated here 
(which is 5K), as consistent with measurements 26]. We also computed the dielectric sus- 
ceptibility, at different temperatures by progressively cooling down the system, from our 
MC simulations via two different approaches: (i) a "direct" method for which the result- 
ing dielectric susceptibility is denoted as ■^ direct and is calculated as the change in polar- 
ization with respect to an applied electric field (with this field practically being oriented 
along the [111] pseudo-cubic direction and having a magnitude of 10 7 V/m); and (ii) the 
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"correlation-function" approaches of Refs. |27l . |28| for which the resulting dielectric sus- 
ceptibility is referred to as x and is provided by the fluctuation-dissipation theorem via 
Xa0 — kjr u a u /3 > ~ < u a >< Up >], where < u a U/3 > denotes the statistical average 
of the product between the a and (3 components of the supercell average of the local mode 
vectors, and where < u a > (respectively, < Up >) is the statistical average of the a- (respec- 
tively, 0-) component of the supercell average of the local mode vectors. N is the number 
of sites in the supercell while V is its volume, ks is the Boltzmann's constant and e D is the 
permittivity of the vacuum. Strikingly, while previous works (see, e.g., Ref. [281 ]) found that 
these two different methods provide nearly identical dielectric susceptibilities in typical fer- 
roelectrics, Fig. 1(a) reveals that it is not the case for disordered BZT: x exhibits a peak 
around T/ ~ 130K, while x ect increases when decreasing the temperature down to Tf and 
then saturates to a plateau for lower temperature. Both the temperature behavior of \ CF 
and the temperature at which \ is maximum are fully consistent with the dielectric ex- 
periments of Ref. 26] in Ba(Zr^Tii_ x )03 relaxors under ac electric fields having frequencies 
ranging between 100Hz and 100kHz. Moreover, the depicted behavior of y direct is exactly the 



one expected for the perfectly static dielectric response of relaxors [23|, lMJ, which allows us 



to identify Tf as the so-called freezing temperature [?-llO| (a freezing temperature ranging 
between 100 and 140K has been reported for BZT systems 29[, in good agreement with 
our value of ~ 130K). Our ■y/ llvect thus provides the static (dc) dielectric response while our 
simulated x corresponds to observed low-frequency dielectric responses of BZT relaxors 
- which is reminiscent of the fact that the susceptibility given by the fluctuation-dissipation 
theorem is nearly the ac-susceptibility in the Edwards- Anderson model of spin glasses 30]. 

It is also important to recall that, while x possesses a peak at Tf, our MC simulations 
indicate that Ba(Zr .5Ti .5)O3 bulk remains macroscopically cubic and non-polar for any 
temperature - as consistent with what is expected for relaxors [1]. Moreover, the temperature 
behaviors of \ CF an d x direct allow the introduction of four different regions, namely (1) 
Region I that concerns temperatures, T, above Th ~ 450 K and for which \ CF and ^/ direct 
can be nicely fitted by the Curie- Weiss formula 221 ] . i.e. they are both directly proportional 
to 1/|T — T \ (where T is practically found here to be very close to — 120K); (2) Region 
II that extends between T* ~ 240 K and Tf, for which x CF increases as the temperature 
decreases but does not follow anymore the Curie- Weiss law, unlike ^/ dtrect ^ (3) Region III 
that is located in-between Tf and T* for which neither x nor x dtrect obey the Curie- 
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Weiss law; and Region IV that occurs for temperature lower than T/, and for which x CF 
decreases as T is reduced while ■y/ lirect is nearly constant there. T;, can be assigned to be 
the Burns temperature while T* can be thought as being the novel critical temperature 



recently found in relaxors 
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12] . The facts that \ CF follows the Curie- Weiss law only for 



temperatures above the Burns temperature and that this Burns temperature is of the order 
of 450K have both been observed in Ba(Zro. 5^0.5)03 Similarly, a negative To Curie 

temperature has also been experimentally extracted in BZT samples 



26j. 



Figure 1 
parameter 



reports the temperature evolution of the so-called Edwards-Anderson (E-A) 



251 ]. QeAi that is calculated as Qea — \(^* u i)t ) ■■> where the inner averaging is 



made over the t Monte- Carlo sweeps while the outer averaging is made over the i lattice 
sites. The behavior of the simulated Qea of BZT bulk versus temperature bears some strik- 
ing resemblance with those predicted by the spherical random-bond-random-field model and 
measured from nuclear magnetic resonance for PbMg!/ 3 Nb2/303 relaxor [if]]. For instance, 
(1) it is small and nearly linearly increases when decreasing the temperature at large temper- 
ature (in Region I); (2) it is large and also linearly increases when decreasing the temperature 
for small temperature (in Region IV); and (3) the qEA-versus-T function is curved upward 
in-between (in Regions II and III). Figure 1(b) also reveals that the temperature behavior 
and values of the overall Edward- Anderson parameter (for any temperature) almost entirely 
originate from the electric dipoles centered on Ti ions. Consequently, the contribution of 
the dipoles belonging to BaZr03 unit cells on the total Edwards- Anderson parameter nearly 
vanishes. Other dramatic differences between local properties associated with Zr versus Ti 
atoms are reported in Fig. 1(c), which shows that not only the average magnitude of the 
local dipoles centered on Zr ions is much smaller than those centered on Ti ions, but its 
temperature behavior is also strikingly different: the dipoles belonging to BaZrOs unit cells 
continuously shrink in average as the temperature is reduced, while the dipoles located inside 
BaTi03 cells suddenly become enlarged when decreasing the temperature below T*. Electric 
diffraction measurements 31] and a model emphasizing the importance of the BaTi03 soft 
mode on the relaxor behavior of BZT {32] are also consistent with our prediction that the 
Ti sites carry much larger dipoles than Zr sites. Moreover, the results from Fig. 1(c) imply 
that, at the lowest temperatures, the Ti atoms displace in average by about 0.16 A, while 
the Zr atoms move by 0.03 A from their cubic, equilibrium positions. Such numbers are 
in remarkable agreement with the values of 0.17 and 0.03 A, respectively, obtained by the 
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first-principles calculations of Ref. [33] for a BZT supercell containing 135 atoms 34]. 

Let us now focus on Figures 2, that display dipolar snapshots within a given (y,z) plane 
at different temperatures, in order to gain a microscopic understanding of relaxors. Figure 
2a reveals that Region I consists of randomly oriented dipoles that are centered on Ti ions 
and that are surrounded by much smaller dipoles located inside BaZrOs cells. As indicated 
by Fig. 2b, some of these Ti sites act as nuclei to the formation of small clusters inside which 
the dipoles begin to be parallel to each other in Region II. We numerically found that the 
polarization of these small clusters in Region II does not automatically lie along a < 111 > 
direction. For instance, the average direction of the local modes inside the bottom cluster 
of Fig. 2b is along an orthorhombic-like direction, namely it is close to [Oil]. Interestingly, 
some of these clusters do not even possess a polarization being parallel to a high-symmetry 
direction in Region II, such as the top cluster of Fig. 2b for which the vector resulting from 
the average of the local modes is equal to (-0.012,-0.052,-0.021) in the (x,y,z) basis - that is a 
triclinic direction. It is interesting to realize that thermal strain measurements 26[ strongly 
suggest that polar nanoclusters can exist in BZT up to ~ 440K, which is consistent with 
our finding of small polar clusters in Region II (that extends up to T& ~ 450 K). 

As the system enters Region III, two novel features occur as it can been inferred from 
Fig. 2c. First of all, more (small) polar clusters form as the temperature is decreased, which 
makes the average magnitude of the Ti dipoles increasing (see Fig. lc). Secondly, some 
of these clusters now possesses a polarization being close to a < 111 > direction, such as 
the left and right clusters displayed in Fig. 2c for which the average local mode is equal 
to (0.043,-0.048,0.043) and (0.034,0.037,0.045), respectively. Note that, while the clusters 
are always formed at Ti sites, they do not necessarily stay at the same sites for different 
temperatures, or even for different MC sweeps at the same temperature, in Regions II and 
III. In that sense, they can be thought of being of dynamical nature rather than being static. 

Below Tf, some of these clusters have considerably grown in size, like the one located 
at the bottom right corner in Figs. 2 d-f. Novel clusters can still form when decreasing 
the temperature in Region IV, such as the one near the bottom left corner of Fig. 2f at 
10K. On the other hand, other clusters are frozen in the sense that they are always located 
at the same region of space and have a polarization that lies along the same direction, 
independently of the temperature and MC sweep in Region IV (see the central and bottom 
right clusters in Figs. 2d-f). While the different clusters possess different numbers of Ti sites 
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and have different overall shapes, they share a common feature in Region IV: they all have a 
polarization being close to one of the eight equivalent < 111 > directions, as consistent with 
the experimental finding that Raman spectra indicate a rhombohedral structure for the polar 
regions at liquid nitrogen temperature in BZT relaxors 35j. As the temperature is reduced 
in Region IV, the matrix possesses Zr-centered dipoles that are significantly shrinking in 
magnitude. This matrix in Regions II, III and IV also possesses individual Ti dipoles that 
are oriented along many different directions, as in Region I. 

To have further insight into the relaxor behavior, let us denote as k max . the vector of the 
first Brillouin zone possessing the largest magnitude of the Fourier transform of the local 
dipoles configuration [361 ] . k max is numerically found to be slightly dependent on the choice 
of the used supercell, but is always a non-highly symmetric vector that is close to neither 
the center nor the boundary of the cubic first Brillouin zone. For instance, in case of a 
12 x 12 x 12 supercell, k max = gf^-(— y + z), where a iat is the lattice constant of the 5- 
atom primitive cell and where y and z are unit vectors along the y- and z-axis, respectively. 
Figure 1(d) reports the temperature evolution of the square of the Fourier transform of 
the local dipoles configuration at k maa; . One can clearly see that, in Regions I and II, this 
quantity is nearly zero. On the other hand, it increases when the temperature decreases 
below T* while still remaining fairly small (around 1.5% of the total spectra gathering 
the Fourier transforms at all possible k-points, at 5K). We interpret such latter results as 
indicative that the different nanopolar regions slightly interact in Regrons III and IV in an 
antiferroelectric-like (or incommensurate [37]- or dipolar-wave-like) fashion. Interestingly, 
antiferroelectricity has been previously reported in some relaxor systems 

Let us now compute the correlation between Ti dipoles (we decided to focus on Ti-Ti 
dipolar correlations because Figs. 2 revealed that the polar clusters only contain Ti sites and 
because Fig. lb shows that the overall Edwards- Anderson parameter mostly only originates 
from Ti dipoles). This correlation is practically defined by 6{r) = — J2i \u-'\\u^+ | > wnere 
the index i runs over all the Nt% Ti-sites of the system and where Ui and Ui +r are the 
local modes in cell i and in the cell centered on the Ti atom (if any) distant from r from 
the cell i, respectively [38|]. A value of 1 (respectively, -1) for 8(r) for a given r would 
indicate that Ti dipoles and their neighboring Ti dipoles distant from r are aligned along 
the same (respectively, opposite) direction. Figure 1(e) shows the value of #(r) for various 
representative r vectors, as a function of temperature. One can see that, in Region I and 
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in average, the Ti dipoles are only (and slightly) correlated with the Ti dipoles centered at 
their first nearest neighboring cells. Such correlation further increases in strength when the 
polar nanoclusters form and become bigger in size and in polarization, as the temperature is 
reduced in Regions II, III and IV. Second and third-nearest neighbors also begin to be more 
correlated in average as the temperature decreases in Regions III and IV. Interestingly, a 
significant anticorrelation (see the negative sign of the correlation) between Ti dipoles being 
distant by 3 lattice constant along the z- (or x- or y-) axis also strongly develops in Regions 
III and IV, which reinforces the previous finding that antiferroelectric-like interactions exist 
within the BZT relaxor system. Note that the Supplemental Material also provides and 
discusses the 6(r) function for all the r-vectors lying in the (y,z) plane at 10K. 

A particularly important feature of our scheme is that we can switch on and off some 
interactions in order to determine their effect on physical properties. We numerically found 
that turning off random fields and random strains does not significantly affect the results 
shown in Fi gs 1 -2, which contrasts with a common belief on the microscopic origins of 



relaxors 



16 



171 1 while being more consistent with models proposed for the homovalent 

n n 

(K,Li)Ta03 relaxor [391 . I40j |. On the other hand, our computations reveal that it is the 
difference in polarizability between Ti and Zr ions that leads to the relaxor behavior in 
BZT. As a matter of fact, annihilating such differences in the simulations leads to (1) 
^direct anc j yCF being equal to each other and continuously decreasing as the temperature 
decreases down to OK (with the system remaining cubic and non-polar); (2) the Edwards- 
Anderson parameter being around ten times smaller than the one depicted in Fig. lb at low 
temperature, and (3) the polar nanoclusters disappearing. It should also be emphasized that 
our simulations results depicted in Figs 1-2 imply that relaxor behavior can occur in BZT 
even if no large chemically-ordered region exists in that system (since our computations were 
performed on disordered solid solutions). Such finding seems to contrast with models recently 
proposed to explain the relaxor behavior of heterovalent Pb(Sc,Nb)03 and Pb(Mg,Nb)03 
alloys 13J], while agreeing with a study 3jj downplaying the role of chemical short-range 
ordering on the formation of polar nanoregions in BZT. In fact, our simulations indicate that 
the relaxor behavior already occurs in disordered BZT solid solutions because some regions 
of space can be more Ti-rich than others because of the random process of assigning sites in 
a disordered solid solution. Such feature bears resemblance with the Anderson localization 
phenomenon for which electronic wave-functions become localized in a region of space (of an 
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overall disordered (A', A") solid solution) being much richer in A' than in A" 41] . Finally, 
we also increased the antiferroelectricity-like interactions (by playing with the so-called J5 
short-range coefficient 



42]). We found that such increase leads to a shift towards higher 
temperature of the peak of x \ i n addition to enhance at low temperature (i) the Edwards- 
Anderson parameter, (ii) the average magnitude of the local modes centered on Ti ions, (hi) 
the square of the Fourier transform of the local dipoles configuration at k maa; and (iv) the 
strength of the anticorrelation between Ti dipoles being distant by 3 lattice constants along 
the z- (or x- or y-) axis. Such findings emphasize the importance of the antiferroelectricity- 
like interactions between Ti-rich nanopolar clusters for the relaxor behavior. 

We therefore hope that our study helps in better understanding the fascinating relaxor 
ferroelectrics. In order to further enhance such understanding, future studies could examine 



the influence of static and dynamic (GHz-THz) electric fields 20j, |21( on the behaviors of 
BZT materials, and determine if the results found here also hold for heterovalent relaxors 
(such as Pb(Sc,Nb)0 3 and Pb(Mg,Nb)0 3 ). 
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FIGURE CAPTIONS 



Figure 1: Temperature dependency of some properties in disordered Ba(Zr .5Ti .5)O3 
solid solutions. Panel (a) shows the average between the three diagonal elements of the 
dielectric susceptibility, as computed from a direct approach (^* rect ; triangles) and from the 
fluctuation-dissipation theorem (x \ dots). Panel (b) displays the total Edwards- Anderson 
parameter, as well as its contributions from cells centered on Ti and Zr ions. Panel (c) reveals 
the magnitude of the local modes centered on Ti and Zr ions. Panel (d) represents the square 
of the Fourier transform of the local modes' configurations at k max . Panel (e) provides the 
6(r) correlation between Ti dipoles for r= ai at z (first nearest neighbor), ai at (y + z) (second 
nearest neighbor), a^ at (x + y + z) (third nearest neighbor), 2ai at z and 3ai at z. The solid line 
in Panel (a) represents the dielectric susceptibility arising from the fit of % (between 500 



and 800K) by the Curie- Weiss law [22]. 



Figure 2: Snapshots of the dipolar configurations in a given (y,z) plane for different 
temperatures. Panels (a), (b), (c), (d), (e) and (f) correspond to temperature of 550K 
(Region I), 250K (Region II), 150K (Region III), 100K (Region IV), 50K (Region IV) and 
10K (Region IV), respectively. Blue colors and red colors indicate that the corresponding 
local modes are centered on Ti and Zr ions, respectively. 
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